################################################################################################
###Replication Code for "Because He's Gay": How Race, Gender, and Sexuality Shape Perceptions of Judicial Fairness
###Journal of Politics
###Authors: Ana Bracic, Mackenzie Israel-Trummel, Tyler Johnson, and Kathleen Tipler
################################


################################
##Instructions:
##This code file is to analyze Republican respondents.
##First run the code file 1_Conjoint Set Up to generate the stacked data file for analysis.
################################



##Restrict to Republicans. Load stacked data. 
cd <- read.csv(file = "stackedconjointdata.csv")
length(cd$id) #40420 observations total

table(cd$Rpid)
reps <- subset(cd, cd$Rpid > 4)
table(reps$Rpid)
names(reps)
length(reps$id) #15100 Republican observations 

reps$race <- as.factor(reps$race)
reps$race <- relevel(reps$race, ref="White")


reps$sexuality <- as.factor(reps$sexuality)
reps$sexuality <- relevel(reps$sexuality, ref="Straight")


reps$lawschool <- as.factor(reps$lawschool)
reps$lawschool <- factor(reps$lawschool, levels = c("Top 10, Public", "Top 10, Private", "#11-25, Private", "#11-25, Public", "#26-50, Private" ,   "#26-50, Public",  "#51-100, Private",   "#51-100, Public", "#101-200, Private",  "#101-200, Public"))

reps$occupation <- as.factor(reps$occupation)
reps$occupation <- relevel(reps$occupation, ref="Law professor")

reps$Rrace <- as.factor(reps$Rrace)
reps$Rrace <- relevel(reps$Rrace, ref="White")





## Load necessary packages for importing the function for robust clustered errors
library(RCurl)

## Import the function from repository
url_robust <- "https://raw.githubusercontent.com/IsidoreBeautrelet/economictheoryblog/master/robust_summary.R"
eval(parse(text = getURL(url_robust, ssl.verifypeer = FALSE)),
     envir=.GlobalEnv)


## New function: summary(modelname, cluster = c("id)) provides you with robust clustered error adjustment


############################################################
### Regressions and Effects Figures
#############################################################


########################################################################
## CREATE FIGURE 2B
########################################################################

## Conservatism of judge
modC.R <- lm(conservative_1 ~ race + gender + sexuality + party + agecat + lawschool + occupation
             + Rwoman + RLGB + Rincome + Rrelig + Rlibcon + Rrace + Redu + Rage + Rknow    
             +factor(vignette_id),
             data = reps)

summary(modC.R)

## use new summary function
summary(modC.R)
modC2.R <- summary(modC.R, cluster = c("id"))



coefsC.R <- c(modC2.R$coef[,1]) #Pull coefficients from model
confintC.R <- c(modC2.R$coef[,2]*1.96) #Create 95 percent confidence interval from standard error 

## Join into data frame
tmpC.R <- data.frame(cbind(coefsC.R,confintC.R))



##Create slimmed data frame to just attributes for plot

tmpC.R2 <- tmpC.R[grep("^race|^genderWoman$|^sexualityGay$|^partyRepublican",row.names(tmpC.R)),]
tmpC.R2$variable <- c("Asian", "Black", "Hispanic", "Native", "Woman", "Gay", "Republican")
tmpC.R2$variable <- as.factor(tmpC.R2$variable)
names(tmpC.R2) <- c("coef", "confint", "variable")


## Load list of variables for plot axis
variablesC.R <- read.csv(file = "conjointvariablesSHORT.csv")


## Input data for plotting
for (i in 1:dim(variablesC.R)[1]){
  if (sum(tmpC.R2$variable == variablesC.R$variable[i], na.rm = T) > 0){
    variablesC.R$coef[i] <- tmpC.R2$coef[tmpC.R2$variable == variablesC.R$variable[i]]
    variablesC.R$confint[i] <- tmpC.R2$confint[tmpC.R2$variable == variablesC.R$variable[i]]
  }
}


## Reorder for plotting
variables2C.R <- variablesC.R

order <- 1:dim(variables2C.R)[1]
variables2C.R <- transform(variables2C.R, variable2 = reorder(variable2, -order) ) 


variables2C.R$labels <- formatC(variables2C.R$coef, format = "f", digits = 3) #to make sure last 0 shows up for one of the labels
variables2C.R$labels[variables2C.R$labels == "0.000"] <- NA


##PLOT

plotC.R <- ggplot(data = variables2C.R, aes(x = coef, y = variable2)) +
  geom_point(aes(colour = group, shape = group), size = 4.5)+
  geom_text(aes(label = labels), nudge_y = 0.4, nudge_x = -0.002, size=6)+
  geom_errorbarh( aes(y = variable2, xmin = coef - confint, xmax = coef + confint, height=.28,colour=group),size=.5)+
  scale_color_manual(values=c("#e76254", "#e76254", "#e76254", "#e76254")) + 
  scale_shape_manual(values=c(16, 16, 16,  16)) + 
  theme_bw() + 
  xlim(-0.1, 0.2) + 
  xlab("Change in Perceived Conservatism")+
  theme(axis.title.y = element_blank()) +
  theme(axis.title.x = element_text(size = 20)) +
  geom_vline(xintercept = 0,size=.5,colour="black",linetype="dotted") +
  theme(axis.text.y = element_text(angle = 0, hjust = 0, color="black", size = 20))  +
  theme(axis.text.x = element_text(size = 20))  +
  theme(legend.position = "none") +
  theme(plot.margin = unit(c(0.1, .5, 0.1, .5), "cm"))


pdf(file = "Fig2B.pdf", height=9, width=8)
plotC.R
dev.off()




########################################################################
## CREATE FIGURE 3B
########################################################################

## Selected to hear case fairly


modS.R <- lm(selected ~ race + gender + sexuality + party + agecat + lawschool + occupation
             + Rwoman + RLGB + Rincome + Rrelig + Rlibcon + Rrace + Redu + Rage + Rknow    
             + factor(vignette_id),
             data = reps)

# use new summary function
summary(modS.R)
modS2.R <- summary(modS.R, cluster = c("id"))


coefsS.R <- c(modS2.R$coef[,1]) #Pull coefficients
confintS.R <- c(modS2.R$coef[,2]*1.96) #Create 95 percent confidence interval from standard error 

## Join into data frame
tmpS.R <- data.frame(cbind(coefsS.R, confintS.R))

## Create slimmed data frame to just attributes for plot
tmpS.R2 <- tmpS.R[grep("^race|^genderWoman$|^sexualityGay$|^partyRepublican",row.names(tmpS.R)),]
tmpS.R2$variable <- c("Asian", "Black", "Hispanic", "Native", "Woman", "Gay", "Republican")
tmpS.R2$variable <- as.factor(tmpS.R2$variable)
names(tmpS.R2) <- c("coef", "confint", "variable")



#input list of variables for plot axis
variablesS.R <- read.csv(file = "conjointvariablesSHORT.csv")

#input data for plotting
for (i in 1:dim(variablesS.R)[1]){
  if (sum(tmpS.R2$variable == variablesS.R$variable[i], na.rm=T)>0){
    variablesS.R$coef[i] <- tmpS.R2$coef[tmpS.R2$variable == variablesS.R$variable[i]]
    variablesS.R$confint[i] <- tmpS.R2$confint[tmpS.R2$variable == variablesS.R$variable[i]]
  }
}

#reorder for plotting
variables2S.R <- variablesS.R


order <- 1:dim(variables2S.R)[1]
variables2S.R <- transform(variables2S.R, variable2 = reorder(variable2, -order) ) 

table(variables2S.R$variable2 , variables2S.R$group)
variables2S.R$labels <- round(variables2S.R$coef, 3)
variables2S.R$labels[variables2S.R$labels == 0] <- NA

#PLOT

plotS.R <- ggplot(data = variables2S.R, aes(x = coef, y = variable2)) +
  geom_point(aes(shape = group, color = group), size = 4.5)+
  geom_text(aes(label = labels), nudge_y = 0.35, size=6)+
  geom_errorbarh( aes(y = variable2, xmin = coef - confint, xmax = coef + confint, height=.28, colour = group),size = .5)+
  scale_color_manual(values = c("#e76254", "#e76254", "#e76254", "#e76254")) + 
  scale_shape_manual(values = c(16, 16, 16,  16)) + 
  theme_bw() + 
  xlab("Change in Judge Preference")+
  xlim(-.17, .15)+
  theme(axis.title.y = element_blank()) +
  theme(axis.title.x = element_text(size = 20)) +
  geom_vline(xintercept = 0, size =.5, colour = "black", linetype = "dotted") +
  theme(axis.text.y = element_text(angle = 0, hjust = 0, color = "black", size = 20))  +
  theme(axis.text.x = element_text(size = 20))  +
  theme(legend.position = "none") +
  theme(plot.margin = unit(c(0.1, .5, 0.1, .5), "cm"))


pdf(file = "Fig3B.pdf", height=9, width=8)
plotS.R
dev.off()


########################################################################
## CREATE FIGURE 4B
########################################################################

## Politically motivated model

modP.R <- lm(political_1 ~ race + gender + sexuality + party + agecat + lawschool + occupation
             + Rwoman + RLGB + Rincome + Rrelig + Rlibcon + Rrace + Redu + Rage + Rknow    
             + factor(vignette_id),
             data = reps)

summary(modP.R)

# use new summary function
summary(modP.R)
modP2.R <- summary(modP.R, cluster = c("id"))


coefsP.R <- c(modP2.R$coef[,1]) #Pull coefficients
confintP.R <- c(modP2.R$coef[,2]*1.96) #Create 95 percent confidence interval from standard error 

## Join into data frame
tmpP.R <- data.frame(cbind(coefsP.R, confintP.R))

## select particular coefficients
tmpP.R2 <- tmpP.R[grep("^race|^genderWoman$|^sexualityGay$|^partyRepublican", row.names(tmpP.R)),]
tmpP.R2$variable <- c("Asian", "Black", "Hispanic", "Native", "Woman", "Gay", "Republican")
tmpP.R2$variable <- as.factor(tmpP.R2$variable)
names(tmpP.R2) <- c("coef", "confint", "variable")



#input list of variables for plot axis
variablesP.R <- read.csv(file = "conjointvariablesSHORT.csv")


#input data for plotting
for (i in 1:dim(variablesP.R)[1]){
  if (sum(tmpP.R2$variable == variablesP.R$variable[i], na.rm=T)>0){
    variablesP.R$coef[i] <- tmpP.R2$coef[tmpP.R2$variable == variablesP.R$variable[i]]
    variablesP.R$confint[i] <- tmpP.R2$confint[tmpP.R2$variable == variablesP.R$variable[i]]
  }
}

#reorder for plotting
variables2P.R <- variablesP.R

order <- 1:dim(variables2P.R)[1]
variables2P.R <- transform(variables2P.R, variable2 = reorder(variable2, -order) ) 

variables2P.R$labels <- formatC(variables2P.R$coef, format = "f", digits = 3)
variables2P.R$labels[variables2P.R$labels == "0.000"] <- NA



#PLOT

plotP.R <- ggplot(data = variables2P.R, aes(x = coef, y = variable2)) +
  geom_point(aes(colour = group, shape = group), size = 4.5)+
  geom_text(aes(label = labels), nudge_y = 0.35, size=6)+
  geom_errorbarh( aes(y = variable2, xmin = coef - confint, xmax = coef + confint, height=.28, colour=group), size=.5)+
  scale_color_manual(values=c("#e76254", "#e76254", "#e76254", "#e76254")) + 
  scale_shape_manual(values=c(16, 16, 16,  16)) + 
  theme_bw() + 
  xlab("Change in Perceived Political Motivation")+
  xlim(-0.06, 0.06)+
  theme(axis.title.y = element_blank()) +
  theme(axis.title.x = element_text(size = 20)) +
  geom_vline(xintercept = 0, size=.5, colour="black", linetype="dotted") +
  theme(axis.text.y = element_text(angle = 0, hjust = 0, color="black", size = 20))  +
  theme(axis.text.x = element_text(size = 20))  +
  theme(legend.position = "none") +
  theme(plot.margin = unit(c(0.1, .5, 0.1, .5), "cm"))


pdf(file = "Fig4B.pdf", height = 9, width = 8)
plotP.R
dev.off()



########################################################################
## CREATE FIGURE 5B
########################################################################



## Trust in supreme court model
modT.R <- lm(trustcourt_1 ~ race + gender + sexuality + party + agecat + lawschool + occupation
             + Rwoman + RLGB + Rincome + Rrelig + Rlibcon + Rrace + Redu + Rage + Rknow    
             +factor(vignette_id),
             data = reps)

summary(modT.R)

# use new summary function
modT2.R <- summary(modT.R, cluster = c("id"))

coefsT.R <- c(modT2.R$coef[,1])
confintT.R <- c(modT2.R$coef[,2]*1.96)

tmpT.R <- data.frame(cbind(coefsT.R, confintT.R))

##Create slimmed data frame to just attributes for plot

tmpT.R2 <- tmpT.R[grep("^race|^genderWoman$|^sexualityGay$|^partyRepublican", row.names(tmpT.R)),]
tmpT.R2$variable <- c("Asian", "Black", "Hispanic", "Native", "Woman", "Gay", "Republican")
tmpT.R2$variable <- as.factor(tmpT.R2$variable)
names(tmpT.R2) <- c("coef", "confint", "variable")


#input list of variables for plot axis
variablesT.R <- read.csv(file = "conjointvariablesSHORT.csv")


#input data for plotting
for (i in 1:dim(variablesT.R)[1]){
  if (sum(tmpT.R2$variable == variablesT.R$variable[i], na.rm = T)>0){
    variablesT.R$coef[i] <- tmpT.R2$coef[tmpT.R2$variable == variablesT.R$variable[i]]
    variablesT.R$confint[i] <- tmpT.R2$confint[tmpT.R2$variable == variablesT.R$variable[i]]
  }
}


#reorder for plotting

variables2T.R <- variablesT.R

order <- 1:dim(variables2T.R)[1]
variables2T.R <- transform(variables2T.R, variable2 = reorder(variable2, -order) ) 

variables2T.R$labels <- round(variables2T.R$coef, 3)
variables2T.R$labels[variables2T.R$labels == 0] <- NA

#PLOT

plotT.R <- ggplot(data = variables2T.R, aes(x = coef, y = variable2)) +
  geom_point(aes(colour = group, shape = group), size = 4.5)+
  geom_text(aes(label = labels), nudge_y = 0.35, size = 6) +
  geom_errorbarh(aes(y = variable2, xmin = coef - confint, xmax = coef + confint, height = .28,colour = group),size = .5)+
  scale_color_manual(values = c("#e76254", "#e76254", "#e76254", "#e76254")) + 
  scale_shape_manual(values = c(16, 16, 16,  16)) + 
  theme_bw() + 
  xlab("Change in Trust") +
  xlim(-.1, .1) + 
  theme(axis.title.y = element_blank()) +
  theme(axis.title.x = element_text(size = 20)) +
  geom_vline(xintercept = 0, size = .5, colour = "black", linetype = "dotted") +
  theme(axis.text.y = element_text(angle = 0, hjust = 0, color = "black", size = 20))  +
  theme(axis.text.x = element_text(size = 20))  +
  theme(legend.position = "none") +
  theme(plot.margin = unit(c(0.1, .5, 0.1, .5), "cm"))


pdf(file = "Fig5B.pdf", height = 9, width = 8)
plotT.R
dev.off()


